function outFeatures = genLBPFs(image)
%get 24 LBP features
LbpImage = image;
waveImage = image;
LBPHist = LBP_riu243(LbpImage);
LBPHist = LBPHist(1:26);
LBPHist = 100*LBPHist/sum(LBPHist);
% LBPHist
% size(LBPHist)
% plot(LBPHist)
firstindex = 5;
lastindex = 20;
% firstindex 
% % lastindex
% if firstindex>10
%     firstindex = 1;
% end
LBPTmp = LBPHist(firstindex:lastindex);
% LBPTmp = LBPTmp./sum(LBPTmp);
% size(LBPTmp)
outFeatures(1) = mean(LBPTmp);
outFeatures(2) = std(LBPTmp);
outFeatures(3) = skewness(LBPTmp);
outFeatures(4) = kurtosis(LBPTmp);
sumx = 0;
sumxx = 0;
% LBPHist = LBPHist./sum(LBPHist);
for j = 1:size(LBPHist,2)
    sumx = sumx+LBPHist(j)*j;
    sumxx = sumxx + LBPHist(j)*j*j;
end
outFeatures(5) = sumx/sum(LBPHist);
outFeatures(6) = sumxx/sum(LBPHist) - outFeatures(5)*outFeatures(5);

%second 6 features
LBPTmp2 = LBPTmp.*LBPTmp;
outFeatures(7) = mean(LBPTmp2);
outFeatures(8) = std(LBPTmp2);
outFeatures(9) = skewness(LBPTmp2);
outFeatures(10) = kurtosis(LBPTmp2);
sumx = 0;
sumxx = 0;
%   Y(find(Y<1500)) = 0;
LBPHist2 = LBPHist.*LBPHist;
for j = 1:size(LBPHist2,2)
    sumx = sumx+LBPHist2(j)*j;
    sumxx = sumxx + LBPHist2(j)*j*j;
end
outFeatures(11) = sumx/sum(LBPHist2);
outFeatures(12) = sumxx/sum(LBPHist2) - outFeatures(11)*outFeatures(11);

[c,s] = wavedec2(waveImage,3,'sym4');
a1 = appcoef2(c,s,'sym4',3);
hh = fspecial('gaussian',3);
a1 = imfilter(a1,hh);
%     h1 = upcoef2('h',H1,'sym4',0);
%     v1 = upcoef2('v',V1,'sym4',0);
%     d1 = upcoef2('d',D1,'sym4',0);
LBPHist = LBP_riu243(a1);
LBPHist = LBPHist(1:26);
LBPHist = 100*LBPHist/sum(LBPHist);
% LBPHist
% size(LBPHist)
% plot(LBPHist)
firstindex = 5;
lastindex = 20;
% firstindex 
% % lastindex
% if firstindex>10
%     firstindex = 1;
% end
LBPTmp = LBPHist(firstindex:lastindex);
% LBPTmp = LBPTmp./sum(LBPTmp);
% size(LBPTmp)
outFeatures(13) = mean(LBPTmp);
outFeatures(14) = std(LBPTmp);
outFeatures(15) = skewness(LBPTmp);
outFeatures(16) = Kurtosis(LBPTmp);
sumx = 0;
sumxx = 0;
% LBPHist = LBPHist./sum(LBPHist);
for j = 1:size(LBPHist,2)
    sumx = sumx+LBPHist(j)*j;
    sumxx = sumxx + LBPHist(j)*j*j;
end
outFeatures(17) = sumx/sum(LBPHist);
outFeatures(18) = sumxx/sum(LBPHist) - outFeatures(17)*outFeatures(17);

%second 6 features
LBPTmp2 = LBPTmp.*LBPTmp;
outFeatures(19) = mean(LBPTmp2);
outFeatures(20) = std(LBPTmp2);
outFeatures(21) = skewness(LBPTmp2);
outFeatures(22) = Kurtosis(LBPTmp2);
sumx = 0;
sumxx = 0;
%   Y(find(Y<1500)) = 0;
LBPHist2 = LBPHist.*LBPHist;
for j = 1:size(LBPHist2,2)
    sumx = sumx+LBPHist2(j)*j;
    sumxx = sumxx + LBPHist2(j)*j*j;
end
outFeatures(23) = sumx/sum(LBPHist2);
outFeatures(24) = sumxx/sum(LBPHist2) - outFeatures(23)*outFeatures(23);